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ABSTRACT 


This  report  concerns  the  behavior  of  cyclic 
liquid  jets  in  pulsating  bubbles.  The  jets  were 
developed  in  bubbles  which  were  situated  on  a vibrat- 
ing platform  that  was  contained  within  a vessel 
capable  of  sustaining  reduced  pressures.  The  entire 
assembly  was  mounted  on  a variable  frequency  vibra- 
tion table  and  driven  at  displacement  amplitudes 
on  the  order  of  a millimeter.  The  reduced  pressure 
enabled  the  jets  to  be  quite  easily  formed  within 
the  vessel;  the  jets  were  formed  in  bubbles  contain- 
ing both  vapor  and  gas,  thus  resulting  in  bubbles 
exhibiting  jet  behavior  over  several  cycles.  The 
bubbles  were  photographed  with  a high  speed  camera 
and  a detailed  analysis  made  of  the  behavior  of  the  jet 
as  a function  of  time.  It  was  found  that  the  bubble 
surface  could  be  well  represented  by  the  first  four 
even  terms  of  a Legendre  polynomial  expansion.  Further, 
a graphic  three-dimensional  display  of  the  bubble 
surface  generated  by  the  expansion  gives  excellent 
agreement  with  the  observed  behavior. 
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I . HISTORICAL  BACKGROUND 


A . Interest  in  Cavitation 

The  existence  of  cavitation  was  first  reported  by 
Osborne  Reynolds*  in  1894,  who  observed  small  cavities  in 
water  flowing  in  a tube  containing  a local  constriction. 
The  name  "cavitation"  was  given  a year  later  by  Froude^ 
when,  in  an  investigation  of  a British  destroyer  which 

O 

failed  to  meet  its  design  speed,  Thcrnycroft  and  Barnaby 
postulated  that  small  cavities  which  formed  in  the  vicin- 
ity of  the  screws  limited  the  thrust.  Strange  destruction 
of  screws  in  fast  steamboats  was  described  by  Silberad 
in  1912,  and  in  1919  Cook^  presented  the  first  attempt 
at  a theoretical  approach  to  cavitation.  Lord  Raleigh5 
improved  and  completed  this  theory  for  spherical  bubble 
col  lapse. 

Corrections  to  the  Raleigh  collapse  for  surface 
tension,  viscosity  and  liquid  compressibility  in  gas- 
filled  bubbles  have  since  been  made.  However,  each  of 
these  studies  has  assumed  spherical  bubble  collapse. 

In  actual  fact,  bubbles  often  collapse  in  close  prox- 
imity to  solid  boundaries  which  gives  an  asymmetrical 
pressure  field  and  results  in  non-spherical  collapse. 


The  actual  mechanism  for  cavitation  damage  has  only 
recently  been  determined.  It  was  first  proposed  by 
Kornfield  and  Suvarov6  that  a liquid  jet  formed  during 
bubble  collapse  was  a possible  damage  mechanism.  More 
recently,  Naude  and  Ellis7studied  nonhemispherical 
cavities  collapsing  in  contact  with  a solid  boundary, 
and  have  shown  that  pressures  caused  by  the  cavity  wall 
striking  the  boundary  are  higher  than  those  resulting 
from  a compression  of  gases  inside  the  cavity,  and  are 
responsible  for  the  damage. 

Recent  research  efforts  have  been  aimed  at  developing 
a model  for  cavitation  bubbles  which  describes  the  jet- 
collapse  phenomenon  .8*9,10 


B . Existing  Methods  of  Observation 

There  are  several  avenues  for  measuring  cavitation 

activity.11  One  can  record  the  radiated  pressure  waves 

as  a cavity  collapses  or  observe  the  flux  of  light  from 

the  luminescence  which  may  be  emitted  during  the  collapse 

of  a transient  cavity.  The  mechanical  process  of  erosion 

can  be  noted  by  measuring  the  weight  loss  due  to  the  cavi- 

12 

tation  effects  on  a sample.  Weissler  has  measured 
cavitation  events  via  a chemical  reaction  in  which  chlorine 
is  released  from  carbon  tetrachloride  dissolved  in  water. 
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The  best  way  to  observe  cavitation  events  is  by  visual 
means.  This  is  no  easy  task  as  cavities  typically 
collapse  in  the  order  of  50  microseconds.  However, 
cavitation  collapse  is  visually  accompanied  by  an  audi- 
ble "click"  and  rapidly  translating  bubbles  are  easily 
observed.  Attempts  at  photographing  collapsing  cavities 
using  ultra-high  speed  cameras  with  framing  rates  of  a 

13 

few  million  FPS  have  recently  been  accomplished.  Ellis 
was  able  to  photograph  bubbles  in  a hydrostatic  pressure 
field  and  in  a gravity  force  field  near  a boundary. 

This  method  of  photography  requires  exotic  photo- 
graphic equipment  and  highly  refined  timing  techniques. 

The  resulting  photographs  are  strongly  backlit  and  the 
structure  of  the  jet  is  difficult  to  discern. 

Currently  this  is  the  only  avenue  for  direct  observa- 
tion of  this  phenomenon. 


II.  INTRODUCTION 


Cavitation  has  been  the  focus  of  many  research  efforts 
and  becomes  very  significant  in  the  light  of  the  damage 
which  it  causes.  Cavitation  damage  costs  the  U.  S.  Navy 
thousands  of  dollars  each  year  for  the  replacement  of 
parts  eroded  by  the  phenomenon.  m 
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Recently  it  has  been  determined  that  the  mechanism 
for  the  damage  is  the  formation  of  a jet  which  impinges 
on  the  surface  as  the  bubble  collapses.  It  is  appropriate 
that  the  focus  of  the  project  is  on  these  jets. 

Any  attempt  to  investigate  further  the  nature  of 
this  jet  and  its  formation  must  naturally  utilize  a 
direct  visual  observation  of  the  bubble  as  it  collapses. 
Let  us  first  examine  the  conditions  under  which  cavitation 
occurs.  Much  of  damage  observed  occurs  in  fluid  machinery 
and  is  caused  by  hydrodynamic  cavitation.  This  form  of 
cavitation  is  due  to  fluid  shear  in  rapid  flow  over 
surfaces.  This  cavitation  phenomenon  is  very  difficult 
to  observe  as  the  bubbles  collapse  in  the  rapidly  moving 
1 i q u i d . 

An  alternative  to  hydrodynamic  cavitation  is  acoustic 
cavitation,  which  is  induced  on  rapidly  vibrating  surfaces 
Cavitation  bubbles  are  formed  and  collapse  in  a relatively 
confined  space  on  the  face  of  the  vibrating  boundary.  Thu 
acoustic  cavitation  is  the  logical  choice  in  attempts  to 
observe  cavity  collapse. 

Cavitation  bubble  collapse  under  normal  conditions 
is  a very  rapid  and  violent  phenomenon.  However,  the 
formation  of  cavities  in  a fluid  can  be  encouraged  by 
increasing  the  vapor  pressure  of  the  liquid  or  decreasing 
the  ambient  pressure.  Thus,  under  these  conditions, 
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cavitation  can  be  induced  at  a lower  level  of  acoustic 
1 4 

pressure.  With  this  fact  in  mind,  the  experimental 
apparatus  described  in  the  following  section  was  designed 
and  built  for  the  purpose  of  creating  cyclic  cavities 
with  only  a small  acoustic  pressure,  whi ch  could  be 
observed  with  relative  ease. 

It  is  the  purpose  of  this  paper  to  first  photograph 
these  bubbles  and  second,  to  develop  a clear,  detailed 
description  of  cavity  collapse. 


III. 
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APPARATUS 


The  feasibility  of  generating  acoustic  cavitation 
at  a reduced  ambient  pressure  for  the  purpose  of  observ- 
ing jet  behavior  in  collapsing  cavities  was  demonstrated 
with  a hastily  constructed  experimental  setup  hereafter 
referred  to  as  the  original  configuration.  After  care- 
ful consideration,  a new  experimental  arrangement  was 
selected,  assembled  and  tested.  It  was  this  system 
. which  was  used  to  generate  the  majority  of  cavities  photo- 

I graphed  during  the  course  of  this  project.  Ultimately 

the  original  configuration  of  equipment  was  modified  to 
improve  its  characteristics  and  was  also  used  to  generate 
some  of  the  cavities  photographed. 
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A.  Original  Configuration  of  Equipment 


1 2 


The  original  configuration  of  the  equipment  is  shown 
in  figure  1 . Acoustic  pressure  is  generated  inside  the 

containment  vessel  by  the  means  of  a paint  shaker  to  which 
it  is  secured.  The  paint  shaker  operates  at  60  Hz  and 
the  amplitude  of  the  vibration  can  be  controlled  by  means 
of  a potentiometer  on  the  front  of  the  unit.  The  cylin- 
der which  holds  the  host  fluid  is  of  aluminum,  14  cm 
in  diameter  and  is  drilled  and  fitted  with  access  nozzles 
to  provide  a means  for  insertion  and  removal  of  the  host 
fluid,  and  a means  for  reducing  the  ambient  pressure 
inside  the  vessel.  A small  platform  is  positioned 
inside  the  cylinder  in  a horizontal  position  to  provide 
a site  for  the  formation  of  cavities  which  would  be 
readily  visible. 

The  amplitude  of  vibration  is  measured  via  a small 
accelerometer  which  outputs  to  a voltmeter.  A linear 
relationship  exists  between  the  output  voltage  of  the 
accelerometer  and  the  amplitude  of  the  vibration  of  the 
system.  When  the  output  from  the  accelerometer  was 
displayed  on  an  oscilloscope,  the  signal  appeared  to 
contain  some  other  frequencies  aside  from  the  60  Hz 
driving  frequency  of  the  shaker  table.  In  order  to 
obtain  a better  idea  of  exactly  how  the  system  was 
vibrating,  the  output  from  the  accelerometer  was  fed  to 


a frequency  spectrum  analyzer  and  the  resultant  spectrum 
displayed  via  a strip  chart  recorder.  It  can  be  seen 
from  figure  2 that  the  largest  component  of  the  signal 
is  the  60  Hz  driving  frequency.  However,  the  first  har- 
monic was  down  only  4 db  and  the  sixth  harmonic  down  only 
6 db  from  the  fundamental. 


B . System  Improvements 

In  an  attempt  to  improve  this  spectrum  and  to  add 
driving  frequency  flexibility  to  the  system,  an  M.B. 
Electronics  Vibration  Table-Power  Amplifier  Unit  was 
employed  in  place  of  the  paint  shaker.15  The  existing 
vessel,  when  filled  with  the  host  fluid,  proved  to  be 
too  massive  for  the  vibration  table  and  could  not  be 
driven  at  a sufficient  amplitude  to  induce  acoustic 
cavitation.  A new, lighter  design  for  the  containment 
vessel  was  needed  if  the  vibration  table  was  to  be 
empl oyed . 

The  first  attempt  at  a new  design  for  the  contain- 
ment vessel  incorporated  an  aluminum  cylinder  of  smaller 
radius  with  plexiglass  viewports.  This  was  mounted 
firmly  on  an  aluminum  plate  which  was  secured  to  the 
moving  element  of  the  vibrating  table.  Due  to  this 
smaller  radius,  cavities  generated  in  this  system  were 
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formation  of  larger  cavities. 

It  was  likely  that  the  higher  harmonics  observed 
in  the  frequency  spectrum  of  the  original  configuration 
were  due  to  higher  modes  of  oscillation  of  the  system. 
It  was  found  that  these  higher  modes  could  be  damped 
considerably  by  adding  additional  cross  members  to  the 
supports  which  clamped  the  cylinder  to  the  vibration 
table.  A frequency  spectrum  analysis  of  the  original 
conf i guration,  i ncl  udi  ng  the  additional  cross  members, 
is  shown  in  figure  3 . Note  that  the  second  harmonic 
which  was  so  large  before  is  now  9 db  down  from  the 
fundamental . 

Both  systems  were  then  available  to  generate  cavi- 
ties for  the  high  speed  photographic  study  of  the  jet 
behavior  in  a collapsing  bubble. 


C.  Calibration  of  Shaker  Table  Amplitude 


For  each  film  run,  the  amplitude  of  the  vibration 
of  the  system  must  be  known  so  that  the  acoustic  pressure 
inside  the  cylinder  can  be  determined.  This  was  accom- 
plished by  calibrating  the  output  voltage  from  the 
accelerometer.  Since  a linear  relationship  exists 
between  the  output  voltage  and  the  vibration  amplitude, 
the  amplitude  of  the  vibration  was  measured  at  various 
dial  settings  and  plotted  vs  the  accelerometer  output 
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The  cyclic  behavior  of  the  bubble  lends  itself  to 
a photographic  technique  which  utilizes  the  stroboscopic 
effect.  However,  in  order  to  get  real  time  photographs 
from  which  data  describing  the  cavity  profile  in  time 
could  be  obtained,  a high  speed  camera  was  employed. 

The  camera  used  was  a Fastax  model  W 1 63269  which 
has  a maximum  filming  capability  of  about  6000  FPS.  An 
adapter,  which  also  served  as  an  extension  tube,  was 
mounted  on  the  camera  so  that  a 50  mm  Nikor  F/1.4  lens 
could  be  utilized.  The  extension  tube-lens  system 
resulted  in  images  on  the  film  which  were  52%  life  size. 
Additional  extension  tubes  were  available.  However, 
their  use  was  prohibited  by  a lack  of  sufficient  depth 
of  focus. 

Due  to  the  magnification  of  the  image  on  the  film 
by  the  extension  tube,  at  short  ranges,  the  depth  of 
focus  for  the  camera  lens  at  full  aperture  was  extremely 
small.  This  necessitates  the  use  of  a large  f-number 
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and  very  bright  illumination.  It  was  subsequently 
determined  that  f8  was  the  minimum  f-number  which  would 
provide  an  acceptable  depth  of  focus,  which  was  approx- 
imately 5 mm. 

The  amount  of  light  necessary  to  properly  expose 
Kodak  4-x  reversal  film  #7277  was  determined  by  trial 
and  error,  using  a filming  rate  of  approximately  5000 
FPS  which  corresponds  to  a camera  voltage  of  180V. 
Illumination  was  provided  by  an  8"  aluminum  reflector 
with  a piece  of  tracing  paper  taped  over  the  front  to 
provide  for  some  diffusion.  Various  combinations  of 
f-numbers  and  light  source  intensities  (varied  by  chang- 
ing the  distance  between  the  containment  vessel  and  the 
reflector  face)  were  tried  which  resulted  in  a final 
determination  of  the  best  combination  of  parameters. 

After  a comparison  of  the  films  resulting  from 
this  test, the  following  set  of  parameters  was  arrived 
at  for  best  results: 


Camera  Voltage 
Lens  Aperture 


180V  (~5000  FPS) 


Illumination  Distance  10  cm  (500W,  8" 

refractor) 

In  a situation  where  the  depth  of  focus  is  so  small, 
the  camera  focus  becomes  critical.  An  accurate  method 
of  insuring  that  the  image  is  focussed  on  the  film 
plane  was  needed.  After  several  attempts  to  adjust  the 
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hole  was  cut  in  the  film  and  covered  with  translucent 


tape.  This  served  as  a "ground  glass"  on  which  the 
image  could  be  focussed  directly  on  the  film  plane. 
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V.  EXPERIMENTAL  OBSERVATIONS  AND  RESULTS 


A . High  Speed  Photography 

Once  the  equipment  and  the  camera  were  ready,  the 
actual  filming  of  the  cyclic  cavities  was  begun.^  Nearly 
2300  feet  of  16  mm  film  was  exposed,  developed  and 
studied.  Consequently,  many  different  types  of  behavior 
were  observed  in  the  collapsing  cavities. 

However,  two  of  the  specific  types  observed  seemed 
to  be  most  prevalent  and  will  be  treated  here.  In  the 
first  case,  a simple  collapse  sequence  with  jet  forma- 
tion is  observed  toward  the  final  stages  of  collapse. 

The  second  case  involves  the  formation  of  an  inverted 
jet  subsequent  to  the  formation  of  the  impinging  jet. 

Sequences  from  actual  films  appear  in  figures  6 and  7. 
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Figure  7.  Sequence  from  Fastax  movie  film  (cont.) 


B . Interpretation  of  Photographs 


Once  the  high  speed  photographs  were  taken  and  a 
section  of  film  was  selected  for  analysis,  a method  for 
digitizing  the  bubble  profile  was  developed. 

First  the  movie  film  was  enlarged  and  printed  frame 
by  frame  to  the  maximum  size  possible  with  a D2V  Omega 
enl arger . 

These  photographs  were  then  overlayed  with  a grid, 
and  points  on  the  cavity  profile  recorded  in  Cartesian 
coordinates.  It  is  assumed  that  due  to  the  wavelength 
of  the  driving  oscillation,  the  cavity  is  essentially  in 
an  infinite  fluid  and  has  axial  symmetry.  Based  on  this 
assumption,  dota  was  taken  only  from  the  first  quadrant 
and  then  reflected  across  the  axis.  The  height  of  the 
cavity  in  the  first  photograph  of  the  sequence  was  then 
set  equal  to  one,  and  the  rest  of  the  data  for  that 
sequence  normalized.  Since  the  magnification  of  the 
lens  is  known,  the  true  height  of  the  bubble  could  be 
determined  from  a direct  measurement  of  the  film. 

Data  was  taken  from  sections  of  film  recording  the 
two  cases.  Case  one  was  without  the  inverted  jet. 

This  data  appears  in  its  entirety  in  Appendix  A1  , and 
edited  versions  in  figures  8 and  9. 

If  the  film  plane  of  the  camera  is  not  parallel 
to  the  vertical  axis  of  the  bubble,  the  actual  profile 
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of  the  bubble  is  not  recorded  but  rather  the  projection 
of  that  profile  onto  a plane  perpendicular  to  the  camera's 
angle  of  view.  (See  figure  10.)  The  effect  of  this  tilt 
must  be  taken  into  consideration  in  any  further  analysis 
of  this  data  and  will  be  developed  more  completely  in 
the  following  section. 

j 

VI  . MATHEMATICAL  INTERPRETATION  OF  RESULTS 

The  data  taken  from  the  photographs  of  the  bubble 
collapsing  in  time  is  in  itself  a mathematical  description 

of  the  jet  behavior.  However,  in  its  present  form,  for  a 

* 

large  volume  of  data,  it  is  rather  unwieldly.  Therefore, 
a more  accessible  analytic  form  for  the  description  of 
the  bubble  profile  is  desirable. 

Previous  efforts  to  describe  the  non-spherical 
collapse  of  cavities  have  been  made  and  the  approach 
here  will  follow  that  of  Naude  and  Ellis.'7 

It  has  been  shown  by  Naude  and  Ellis  7 that  for  a 
cavity  which  is  symmetric  about  a line  normal  to  a solid 
plane  boundary  and  assuming  inviscid,  incompressible 
fluid,  a constant  pressure  inside  the  cavity,  and  neglect- 
ing surface  tension,  adhesion  and  the  effect  of  gravity, 

i 

that  the  problem  reduces  to: 
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V24>  ( r ,0  , t ) = 0 , 


where  6 is  the  polar  angle  measured  from  the  vertical, 
and  t is  a scalar  potential  function  such  that  V4>  = q, 
the  velocity  vector. 

Solving  this  equation  with  the  appropriate  boundary 
conditions,  Ellis7  has  shown  that 

. *411  ♦ t e2n(t)  p2n(cose) 


where  P2n(cos0)  are  the  Legendre  polynomials  and  q> 2 n ( t ) 
are  the  time  dependent  coefficients  in  the  expansion. 

The  approach  leads  to  an  expression  for  the  bubble 
surface  of  the  form 


oo 

R(e,t)  = a2n(t)P2n(cose) 


Thus  the  problem  is  reduced  to  first  determine  the  best 
set  of  coefficients  which  describe  each  bubble  profile 
and  to  second,  evaluate  these  coefficients  and  determine 
the  validity  of  truncating  the  infinite  series  at  a 
finite  number  of  terms. 

Before  a rigorous  analysis  of  the  data  can  take 
place,  the  question  raised  earlier  of  the  camera  view 
angle  must  be  treated.  In  the  following  discussion  the 
effect  of  a small  rotation  in  camera  angle  of  the  coeffi- 
cients of  the  first  four  even  Legendre  polynomials  will 
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be  considered. 

Refer  once  aqain  to  figure  10.  If  the  film  plane  of  the 
camera  is  tilted  at  an  angle  <p  with  respect  to  the  vertical,  the 
actual  image  on  the  film  will  be  the  projection  of  the  object 
back  onto  a plane  parallel  to  the  film  plane.  This  leads  to  a 
systematic  reduction  in  the  height  of  the  image  at  any  point 
which  is  now  given  by  y(o)  = y(0)cos<{>  . This  gives  rise  to 

a change  in  each  of  the  Legendre  coefficients,  which  are  now 
given  by: 

*0  = ^00  + ^2oH\o+^60 

3 2 " a2+  a0A02+a2^22+a4^42  +36A62 

a4  _ <*4  + ^ o Aq  4+3  2 A2  4 +a,A4  4 +^64 
a6  “ a6+  ^ a0^06+a2^26+a4^46+a6^66 

For  angles  less  than  10°  this  represents  a very  small 
change  in  the  coefficients.  This  is  described  in  detail  in 
Appendix  A2. 
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The  method  of  least  squares  was  used  in  determining  the 
coefficients  for  the  first  four  even  Legendre  polynomials  in 
cos9 , which  best  fit  each  of  the  profiles  analyzed. 

As  developed  earlier  the  bubble  surface  can  be  described  by: 


R(e,t)  = a2n(t)P2n(cos0) 


n=0 


The  coefficients  ^(t)  were  determined  from  the  following 
analysis: 

N 

Let  S(a „,a  .a  ,a J = V"*  {^(6, a ,a  ,a  ,a  )-r.}2 

0’  2 I*  6 / -i  1 0 2 4 6 1 

i=l 


where  N is  the  number  of  points  which  describes  the  bubble  surface. 
Thus  S represents  the  sum  of  the  squares  of  the  differences  between 
the  function  f and  the  actual  surface. 

In  order  to  fit  f to  the  points  which  describe  the  surface. 


S must  be  a minimum  which  implies 


al=  0 . 


Where 


S(a0,a2,a4,a6)  = a0P0(cose)  + a2P2(cos0) 

+ a4P4(cos0)  + a6P6(cos0) 


This  results  in  four  simultaneous  equations  in  the  four  coeffi- 
cients 


a EP7 

0 0 


9o  £P 0^2 
3qTP oP 4 


anEPnP 
0 0 


0 


+ a IP  P + a IP  P + a IP  P = Er.P 

202  404  60S  10 

+ a2EP|  + a4EP2P4  + a6EP2PG  = Er.P2 

+ a2EP2P4  + a4EP*  + a6EP4P6  = Er.P„ 

+ a,EP,P,  + a, EP, Pc  + acEP^  = Er.Pc 

226  446  66  16 


Representing  this  set  of  equations  as  a matrix  equation  we  have 
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These  equations  can  now  be  written  in  matrix  notation: 

S ( 4x4  ) x A( 4x1 ) = C ( 4x1  ) 

Solving  for  A : 

A(4xl ) = S'1 (4x4)  x C ( 4 x 1 ) 

Where  A(4,l)  are  the  coefficients  which  best  fit  the 
polynomials  to  the  data. 

These  matrix  manipulations  lend  the  problem  to 
numerical  solution  by  computer.  A program17  was  written 
to  determine  the  coefficients  and  then  to  generate  a 
plot  of  the  surface  using  those  coefficients.  The  pro- 
gram also  computes  the  x2  as  a test  of  the  accuracy 
of  the  fit. 

This  fit  of  the  Legendre  polynomials  is  meaningful 
only  for  bubble  profiles  in  which  the  jet  has  not  yet 
struck  the  boundary.  At  that  point,  the  surface  is  no 
longer  a single  valued  function  in  r and  e,  and  this 
analysis  is  no  longer  valid.  However,  once  the  coeffi- 
cients are  determined  for  the  first  ten  or  so  steps 
which  describe  the  collapse  up  to  the  frame  just  before 
the  jet  impinges  on  the  surface,  a short  extrapolation 
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two  sets  of  data  describing  the  two  cases.  Once  the 
coefficients  were  solved,  a cubic  spline  fit  was  used 
to  fill  in  additional  coefficients  for  the  times  in 
between  frames.-'0  These  two  resulting  sets  of  coeffi- 
cients appear  graphically  in  figures  11  and  10  where 
the  actual  coefficients  for  each  frame  analyzed  can  be 
found  in  Appendix  5. 

Figures  13  through  15  show  examples  of  the  normal- 
ized data  and  the  resulting  fit.  The  values  for  the  x? 
Goodness -of-fi t Test  for  both  cases  are  illustrated 
in  figurelF  and  appear  in  tabulated  form  in  Appendix  6. 

A sophisticated,  three  dimensional  graphic  display 
system  at  the  Naval  Academy  provided  an  excellent 
opportunity  to  portray  the  history  of  the  bubble  collapse 
via  the  coefficients  which  were  generated  earlier.  A 
photograph  of  the  system  is  shown  in  figure  17. 

A program  was  prepared  for  the  system's  PDP-11  mini- 
computer to  compute  and  display  three  dimensional  data 
for  the  bubble  surface  from  each  set  of  four  coeffi- 
cients.19 As  the  computer  generated  the  data  for  each 
surface  in  real  time,  it  was  displayed  in  sequence 
on  a cathode  ray  tube.  Thus,  the  collapse  of  the  bubble 
was  reconstructed  visually  using  the  coefficients  which 
had  been  generated  earlier.  From  this  picture  screen, 
a movie  was  taken  which  gives  both  a profile  view  and  a 
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Figure  16.  X Goodness-of-fit  Test. 


30°  look  down  angle  for  each  collapse  case. 

A second  program  was  then  written  to  provide  for 

hard  copy  generation  of  the  3 dimensional  data  on  a 

20 

Xynetics  flat  bed  plotter.  Selected  frames  from  the 
collapse  sequence  could  then  be  reproduced  on  paper 
for  more  permanent  display  purposes.  Examples  of  these 
appear  in  figures  18  and  19. 

VII.  SUMMARY  AND  CONCLUSIONS 

It  is  helpful  now  to  review  the  basic  aspects  of 
this  project.  First,  the  equipment  for  the  generation 
of  cyclic  pulsating  cavities  at  a reduced  ambient 
pressure  was  designed,  built,  tested,  and  modified. 

Next,  it  was  demonstrated  that  cyclic  cavities  could 
be  localized  on  the  horizontal  plate  inside  the  contain- 
ment vessel.  The  behavior  of  the  oscillating  bubbles 
was  then  observed  and  recorded  with  a Fastax  movie 
camera.  After  a reliable  set  of  photographic  parameters 
was  determined  through  testing  of  the  camera  on  the 
actual  system,  films  were  taken  for  the  purpose  of 
obtaining  a digitized  description  of  the  bubble  pro- 
file during  collapse.  Following  this,  the  data  was 
analyzed  via  a computer  program  which  fit  the  data  to 
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the  first  four  even  Legendre  polynomials  so  that  an 
analytic  description  of  the  bubble  profile  in  time 
could  be  obtained. 

The  above  outlined  method  of  investigation  gave 
excellent  results.  As  can  be  seen  from  the  x2  Goodness- 
of-fit  Test,  the  fit  was  remarkably  good.  Also,  by  a 
short  extrapolation,  the  coefficients  which  describe 
the  complete  history  of  the  collapse  were  obtained. 

These  results  are  significant  in  that  this  type  of 
analysis  of  raw  data  has  resulted  in  the  generation  of 
Legendre  coefficients  for  the  collapse  in  time.  It  is 
at  this  point  that  future  theoretical  approaches  to 
the  analytical  behavior  of  collapsing  cavities  will  be 
compared  with  theory. 

Finally,  through  the  use  of  the  Naval  Academy's 
PDP-11,  / F-f$  Picture  System,  a 3 dimensional  display 
was  generated  of  the  cavity  collapse  in  time  as  described 
by  the  coefficients  obtained  earlier. 

This  display  provided  for  a visual  comparison  between 
the  actual  photographs  of  the  collapsing  bubble  and  the 
resulting  fit  of  the  Legendre  polynomials. 

In  conclusion,  this  graphic  representation  of 
cavity  collapse,  as  well  as  the  detailed  visual  record 
obtained  from  the  high  speed  photographs,  is  unique 


with  regard  to  the  existing  literature  and  represents 
a significant  advance  in  cavitation  research. 

Some  preliminary  results  of  this  project  were 
reported  at  the  Austin  meeting  of  the  Acoustical 
Society  of  America  and  were  quite  well  received. 
Requests  for  more  details  concerning  the  collapse 
sequence--the  jet  behavior  has  never  been  seen  so 
clearly  before--have  been  received  from  a number  of 
researchers.  Current  plans  are  to  send  a copy  of 
edited  portions  of  the  filmed  jet  collapse  sequence  to 
these  researchers  for  possible  collaboration. 
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0.235294 
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0.552941 

0.352941 

0.552941 

0.411765 

0.547059 

0.470588 

0.535294 

0.529412 

0.523529 

0.588235 
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0.452941 

0.705882 
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0.764706 
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0.8235^9 

0.235294 

0.832353 

0.176471 
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5.33235 
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0.470538 
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0.432353 

0.235294 

0.511755 

0.294113 

0.517647 

0.352941 

0.517547 

0.411765 

0.517647 

0.470533 

0,508324 

0.529412 

0.470533 

9.647059 

0.494113 

0.533235 

0.435294 

0.705832 

0.411765 
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0.329599 

0.823529 

0.235294 

0,332353 
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0.176471 

0.464706 
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0.470583 

0.294113 

0.470583 

0.352941 
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0.411765 

0.464706 

0,470588 

0.447059 

0.529412 

0.408824 

0.647059 

0.435294 
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0.364706 

0.705332 

0.305882 

0.764706 

0.235294 

0.823529 

0.176471 
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0.364706 

5.88235  E 
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i° 

0.332353 

FRAME  8 
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0.352941 
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0.383235 

5.88235  E-2 

0.411765 

0.117647 

0.423529 

0.176471 

0.429412 

0.235294 

0.429412 

0.294118 

0.423529 

0.352941 

0.417647 

0.411765 

0.411765 

0.470583 

0.333295 

0.529412 

0.376471 

0.588235 

0.347059 

0.647059 
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0.705882 

0.235294 

0.764706 

0.176471 

0.794118 

0.117647 

0.817647 

5.33235  E 

:-2  0.323529 

0 

0.323529 
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frame  11 
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5.38235  E-2  0.470533 
2.94113  E-2  0.494113 
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i 

0,2  '7 757 

f) 

0.274113 
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0,21411 • 

0,206382 
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0.714706 
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0 
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0.147059 
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0.176471 

0.329412 

0.235294 
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0.264706 

0.235294 

0.294113 

0.176471 

0.3 

0.141176 

0.288235 

0.114706 

0.264706 

7.64706  E- 
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APPENDIX  A2 

Correction  for  camera  view  angle. 
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r ( e j ) s i n e j = r'(e)sine 

tan9j  = tanecosij) 

r(0x  )cose  j cosif  = r'(e)cos<i> 

r'(0)  = r(e  j )^si  n201  + cos2©  x cos2$ 

r 1 ( 0 ) = rf  tan'Mtanflrns^)  j.  /Ta~n^  e cos  2 j>  + cos2<fr 

\l +tan20cos24>  1 + 1 a n 2qcos2<}) 


First  evaluate 


tan2ecos2(|)  + cos2(f> 
l+tan2ecos2<)>  l+tan20cos2<t> 

cos2i|){l+tan20} 

1 + tan20cos2c(i 

_ cos2j> 

cosz0  + si  n2ecos2<f> 

_ COS2(|)  

cos^esin'^  + cos^ij) 


1 + coszetanz!t> 


r 1 ( 0 ) = rltan"1  (tan  cos4> ) } 
yi  + cos2©tan2<f> 

Let  r(©)  = a + a P (cose)  + a P (cose)  + a P (cose) 

0 2 2 4 4 6 6 

r ( 0 ) = 1 akPk(cose)  where 

k=0,2,4,6 

and  let  r‘(e)  = 1 akPk(cose) 

k=0,2,4,6 


Ea'P  (cos6)  = Ea,  kjcos f tan  ( tanPcostfa)  ,1 

(1  + cos20tan24)) l' 2 
k=0,2 ,4 ,6  k=0,2,4,6 


Multiply  both  sides  by  P .(cose)sinode  and  integrate  from  Ott 

J 


using  the  prooerty 


f P, (cos0)P . (cose )sin0de  = — — 5kj 
J K J 2j+l 


7T 

Eak  ^*Pk(coso)Pj(coso)sinede  = Eak  j"  ?k  cos{tan  1 ( tanecosg>) } 

„ (1  + cos^Otan^^'T7 


k=0,2,4,6 


k=0,2,4 ,6 


x P.(cose)sinnd0 
0 


Also  notice  that 
cosa 


1 


So 


V 


1 + tan2a 


cost  tan~ 1 (tanscostp)}  = 


(1  + tan2{tan ’’ (tan0cos4>) }) 1/2 


(1  + tan20cos2<}>) iy^2 


And 


a ^ ™ 
J 2 


:a  JV I 

n-'o  1 + tan2( 


’o  1 + tan29cos2<j))  */2 
(1  + cos20tan24>) */2 
k=0,2,4 ,6 


X P . (coso)sined9 
J 
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The  integral  becomes 
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So 


aJ 


^(|i) 


1 akAkj 
k=0,2,4,6 


We  now  must  evaluate 
Use  the  relation 

dP. (x)  , , , 

(1-x2)  _ L = kPk_-, (x)  - kxPk(x) 

i 

Ay  * /{ kxPk.,(x)  - kx'P^lx)  - Pk(x)x“)Pj(x)dx 
" I 1 

Akj  =f  kxP|>]P.(x)dx  - (k+l)/x2Pk(x)Pj(x)dx 

Also  since  K,j  = 0,2,4, 6 Pk  is  even  and  P|(_1  is  odd: 


am  - 2 1 k y 
1 0 

xPk-i (x)pj(x)dx  - 

(M>) 

0 

’x2Pk(x)Pj(x)dx J 

A0  0=  ~2/3 

A02  = -4/15 

o 

■r 

II 

0 

o 

<T> 

II 

0 

A20=  8/15 

A22  = -2/21 

A24 

-24/105 

A26 

0 

A4<T  0 

A22  = 32/105 

A24= 

-38/693 

A26  = 

-100/429 

> 

cn 

o 

II 

O 

A = 0 
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40/143 
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-82/2145 
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Program  which,  given  the  normalized  data, 

L I computes  the  coefficients  of  the  first  four 

even  Legendre  polynomials  in  a least  squares 


' 
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TRIDENT 


100  FILE  #3:"COEFG" 

101  FILE#4 : "CHIG" 

110  LET  B=4 

120  DIM  S(4,4),M(4,4),A(4,1),C(4,1),P(4) 

130  PRINT"INPUT  THE  DATA  FILE  AND  AN  EMPTY  FILE" 

140  INPUT  A$,B$ 

150  FILE#1:A$ 

160  FILE#2 : B$ 

170  SCRATCH#2 

180  INPUT# 1 :X,Y 

190  PRINT#2:X;",";Y 

200  LET  L=L+1 

210  LET  R=SQR(X*X+Y*Y) 

220  LET  Q=3. 14159/2 
230  IF  Y=0  THEN  250 
240  LET  Q=ATN(X/Y) 

250  LET  P(1)=1 

260  LET  P(2)=(3*COS(2*Q)+1 )/4 

270  LET  P(3)=(35*COS(4*Q)+20*COS(2*Q)+9)/64 

280  LET  P(4)=(231*C0S(6*Q)+126*C0S(4*Q)+105*C0S(2»Q)+50)/512 

290  FOR  N=1  TO  B 

300  FOR  M=1  TO  B 

310  LET  S(N,M)=P(N)*P(M)+S(N,M) 

320  NEXT  M 

330  NEXT  N 

340  FOR  F=1  TO  B 

350  LET  C(F , 1 )=R*P(F)+C(F, 1 ) 

360  NEXT  F 

370  IF  MORE#1  THEN  180 
380  MAT  M=INV(S) 

390  MAT  A=M»C 
400  MAT  PRINT  A 

410  PRINT#3:A(1,1);",";A(2,1);",";A(3,1);",n;A(4,l) 

420  PRINT#2:1E37;",";1E37 

430  FOR  Q=0  TO  3.14159/2  STEP  3.14159/100 

440  LET  P0=1 

450  LET  P2=(3*COS(2*Q)+1)/4 

460  LET  P4=(35*C0S(4*Q)+20*C0S(2*Q)+9)/64 

470  LET  P6=(231*COS(6»Q)+126*COS(4»Q)+105*COS(2*Q)+50)/512 

480  LET  R=A( 1 , 1 )*PO+A(2, 1 )*P2+A(3» 1 )*P4+A(4, 1 )*P6 

490  IF  H=5  THEN  620 

500  LET  X=R*SIN(Q) 

510  LET  Y=R*COS(Q) 

520  PRINT#2:X;",";Y 

530  NEXT  Q 

540  LET  H=5 

550  RESET# 1 

560  INPUT# 1:X,Y 

570  LET  R1=SQR(X*X>Y*Y) 

580  LET  Q=3. 14159/2 
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TRIDENT  (continued) 

590  IF  Y=0  THEN  6 10 
600  LET  Q=ATN(X/Y) 

610  GO  TO  440 

620  LET  K=K+(R1-R)~2/R1 

630  IF  MORE#1  THEN  560 

640  LET  K=K/L 

650  PRINT"X~2  =";K 

655  PRINT#4:K 

660  END 
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APPENDIX  A4 


Program  which  utilizes  a cubic  spline  fit 
and  computes  additional  points  in  between 
given  points. 


- 
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-i-  L 

SPLINFIT 

100  LIBRARY  "L.IG***:TELSUB","L.IG«”:TYPSUB","L.IG«”:CRTSUB" 

110  LIBRARY  "L . IG*** : COMSUB" 

120  ' THIS  IS  THE  INPUT  PORTION  OF  THE  PROGRAM 
130  PRINT 
140  PRINT 

150  DIM  P(4f101),U(4,101),N(101,4),B(4,101)fL(101),D(101) 

160  PRINT  "IS  THE  DATA  IN  AN  INPUT  FILE  IN  STANDARD  GRAPHICAL  FORMAT"; 

170  INPUT  X$ 

180  IF  X$="NO"  THEN  220 

’ 190  PRINT  "INPUT  FILE  NAME"; 

200  INPUT  B$ 

210  FILE  #2:B$ 

220  PRINT  "ENTER  NAME  OF  OUTPUT  FILE"; 

230  INPUT  D$ 

240  FILE  #1 :D$ 

250  SCRATCH  #1 

260  PRINT  "SPACE  OR  PLANE  CURVE"; 

270  INPUT  B$ 

280  IF  X$="NO"  THEN  560 
290  IF  SEG$(B$, 1 ,3)="PLA"  THEN  500 
300  IF  SEG$(B$, 1 ,3)="SPA"  THEN  330 
310  PRINT  "ILLEGAL  OPTION,  RETYPE"; 

320  GO  TO  270 
330  FOR  1=1  TO  4 
340  FOR  J=1  TO  101 
350  LET  P(I,J)=0 
360  LET  U(I,J)=0 
370  LET  N(J,I)=0 
380  LET  B(I, J)=0 
390  LET  L(J)=0 
400  LET  D(J)=0 
410  NEXT  J 
420  NEXT  I 
430  LET  J=0 

440  IF  SEG$(B$,1,3)="PLA"  THEN  500 
450  LET  J=J+1 
460  IF  END#2  THEN  550 
470  INPUT  #2:P(1,J),P(2,J),P(3*J) 

480  IF  P(1,J)=1E37  THEN  550 
490  GO  TO  450 
500  LET  J=J+1 
510  IF  END  #2  THEN  550 
520  INPUT  #2;P(1,J),P(2,J) 

530  IF  P(1,J)=1E37  THEN  550 
540  GO  TO  500 
550  IF  J=1  THEN  980 
560  IF  SEG$(B$,1,3)="SPA"  THEN  620 
570  IF  SEG$(B$,1,3)="PLA"  THEN  600 
580  PRINT  "ILLEGAL  OPTION,  RETYPE"; 

590  GO  TO  270 
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SPLINFIT  (continued) 

600  LET  S=2 

610  GO  TO  630 

620  LET  S=3 

630  LET  N=J-1 

640  IF  X$="N0"  THEN  660 

650  GO  TO  760 

660  PRINT  "INPUT  NUMBER  OF  POINTS  TO  BE  ENTERED"; 

670  INPUT  N 

680  PRINT  "INPUT  POINT  X,Y,  AND  (IF  NECESSARY),  Z COMPONENTS" 

690  PRINT  "ONE  SET  AT  A TIME" 

700  FOR  J=1  TO  N 
710  IF  S=3  THEN  740 
720  INPUT  P( 1 , J) ,P(2, J) 

730  GO  TO  750 

740  INPUT  P(1 , J) ,P(2,J) ,P(3,J) 

750  NEXT  J 
760  LET  R=R+1 
770  PRINT 
780  PRINT 

790  PRINT  "CURVE";R;":" 

800  LET  Z9=0 
810  PRINT 

820  PRINT" INPUT  INITIAL  END  CONDITION"; 

830  INPUT  A$ 

840  PRINT 

850  IF  SEG$(A$,1,3)=  "CYC"  THEN  930 
860  IF  SEG$(A$,1,3)=  "ANT"  THEN  930 

870  CALL  "ENDPNT":S,N,P(  , ),N(  , ),B(  , ),D(  ),L(  ),Z,U(  , ),A$,F$,Z9 

880  IF  Z9=1  THEN  800 

890  IF  SEG$(A$,1,1)=  "Q"  THEN  930 

900  IF  SEG$(F$, 1 , 1 )=  "Q"  THEN  930 

910  CALL  "MATOPS":S,N,P(  , ),N(  , ),B(  , ),D(  ),L(  ),U(  , ) 

920  GO  TO  940 

930  CALL  "MATINV":S,N,P(  , ),N(  , ),B(  , ),D(  ),L(  ),U(  , ),A$,F$,Z 

940  CALL  "CURGEN":S,N,P(  , ),L(  ),Z,U(  , ),#1 

950  LET  W=Z*(N-1 )+2 

960  IF  X$="NO"  THEN  980 

970  GO  TO  330 

980  CALL  "GRAPH" :W,#1,S,D$ 

990  PRINT  "DO  YOU  WANT  ANOTHER  PLOT  OF  THE  SAME  CURVE"; 

1000  INPUT  A$ 

1010  IF  A$="NO"  THEN  1060 

1020  RESET  #2 

1030  SCRATCH  #1 

1040  LET  R=0 

1050  GO  TO  330 

1060  END 

1070  REM  SETS  UP  MATRIX  END  PIS  FOR  REL,CLA,P,Q  CONDITIONS 

1080  SUB  "ENDPNT":S,N,P(  , ),N(  , ),B(  , ),D(  ),L(  ),Z,U(  , ),A$,B$,Z9 

1090  IF  SEG$(A$, 1 , 1 ) = "P"  THEN  1300 
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SPLINFIT  (continued) 

1100  IF  SEG$(A$, 1 , 1 ) = "Q"  THEN  1320 

1110  IF  SEG$(A$, 1,3)  = "REL"  THEN  1270 

1120  IF  SEG$(A$, 1,3)=  "CLA"  THEN  1160 

1130  PRINT"ILLEGAL  INITIAL  CONDITION" ;A$; "RETYPE" 

1140  LET  Z9=1 
1150  GO  TO  2040 

1160  PRINT  "INPUT  INITIAL  END  TANGENT  VECTOR  COMPONENTS"; 

1170  IF  S=3  THEN  1200 
1180  INPUT  U( 1 , 1 ) , U(2f 1 ) 

1190  GO  TO  1210 

1200  INPUT  U( 1 , 1 ) ,U(2, 1 ) ,U(3, 1 ) 

1210  LET  N( 1 ,2)=1 
1220  LET  N( 1 ,3)=0 
1230  FOR  K=1  TO  S 
1240  LET  B(K, 1 )=U(K, 1 ) 

1250  NEXT  K 
1260  GO  TO  1320 
1270  LET  N( 1 ,2)=1 
1280  LET  N(1,3)=.5 
1290  GO  TO  1320 
1300  LET  N( 1 ,2)=1 
1310  LET  N(1,3)=1 
1320  FOR  J=1  TO  N-1 
1330  IF  S=3  THEN  1370 

1340  LET  D( J)=SQR( (P( 1 , J+1 )-P( 1 , J) )~2+(P(2, J+1 )-P(2,J))“2) 

1350  LET  L(J)=D(J) 

1360  GO  TO  1390 

1370  LET  D( J)=SQR( (P( 1 , J+1 )-P( 1 , J) )~2+(P(2, J+1 )-P(2, J) )~2+(P(3, J+1 )-P(3|J))~2) 
1380  LET  L(J)=D(J) 

1390  NEXT  J 

1400  IF  SEG$ ( A$ , 1 , 3 ) = "CLA"  THEN  1600 
1410  IF  SEG$ ( A$ , 1 , 1 ) = "P"  THEN  1480 
1420  IF  SEG$(A$ ,1,1)=  "Q"  THEN  1520 
1430  CHANGE  A TO  A$ 

1440  FOR  K=1  TO  S 

1450  LET  B(K,1)=(3/(2*L(1)))«(P(K,2)-P(K,1)) 

1460  NEXT  K 
1470  GO  TO  1600 
1480  FOR  K=1  TO  S 

1490  LET  B(K, 1 )=(2/L( 1 ) )*(P(K,2)-P(K, 1 )) 

1500  NEXT  K 
1510  GO  TO  1600 

1520  REM  TREAT  Q-SPLINE  CONDITIONS 
1530  LET  S1=L( 1 )/L(2) 

1540  LET  N( 1 , 1 )=1 
1550  LET  N(  1 f2)  = 1-Sl‘*2 
1560  LET  N(1,3)=-Sr2 
1570  FOR  K=1  TO  S 

1580  LET  B(K, 1 )=(2/L( 1 ) )*(-P(K, 1 )+( 1+S1*‘3)#P(K,2)-(S1“3)#P(K,3)) 

1590  NEXT  K 
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SPLINFIT  (continued) 


1600  PRINT 

1610  PRINT"INPUT  FINAL  END  CONDITION"; 

1620  INPUT  B$ 

1630  IF  SEG$(B$,1,1)="P"  THEN  1870 
1640  IF  SEG$(B$,1,1)="Q"  THEN  1930 
1650  IF  SEG$(B$,1,3)="REL"  THEN  1810 
1660  IF  SEG$(B$,1f3)="CLA"  THEN  1700 
1670  PRINT  "ILLEGAL  OPTION,  RETYPE"; 

1680  GO  TO  1610 
1690  ' TEST  POINT 

1700  PRINT  "INPUT  FINAL  END  TANGENT  VECTOR  COMPONENTS" ; 

1710  IF  S=3  THEN  1740 
1720  INPUT  U( 1 ,N) ,U(2,N) 

1730  GO  TO  1750 

1740  INPUT  U(1,N),U(2,N),U(3,N) 

1750  LET  N(N,1)=0 
1760  LET  N(N,2)=1 
1770  FOR  K=1  TO  S 
1780  LET  B(K,N)=U(K,N) 

1790  NEXT  K 
1800  GO  TO  2000 
1810  LET  N(N, 1 )=2 
1820  LET  N(N,2)=4 
1830  FOR  K=1  TO  S 

1840  LET  B(K,N)=(6/L(N-1 ))*(P(K,N)-P(K,N-1 )) 

1850  NEXT  K 
i860  GO  TO  2000 
1870  LET  N(N, 1 )=1 
1880  LET  N(N,2)=1 
1890  FOR  K=1  TO  S 

1900  LET  B(K,N)=(2/L(N-1))«(P(K,N)-P(K,N-1)) 

1910  NEXT  K 

1920  GO  TO  2000 

1930  LET  S2=L(N-1 )/L(N-2) 

1940  LET  N(N,1)=-S2“2 
1950  LET  N(N,2)=1-S2~2 
I960  LET  N(N,3)=1 
1970  FOR  K=1  TO  S 

1980  LET  B(K,N)=(2/L(N-1))»((S2"3)*P(K,N-2)-(1+S2“3)*P(K,N-1)+P(K,N)) 
1990  NEXT  K 
2000  PRINT 

2010  PRINT"NUMBER  OF  INTERMEDIATE  POINTS  PER  SPAN"; 

2020  INPUT  Z 
2030  LET  Z=Z+1 
2040  SUBEND 

2050*  THE  FOLLOWING  SETS  UP  THE  INTERIOR  PORTION  OF  THE  MATRIX 
2060'  NOTE-  HERE  I USE  A DIAGONAL  MATRIX  TO  SAVE  STORAGE  SPACE,  AND 
2070'  EMPLOY  APPROPRIATE  ALGORITHMS  INSTEAD  OF  MAT  FUNCTIONS. 

2080  SUB  "MATOPS":S,N,P(  , ),N(  , ),B(  , ),D(  ),L(  ),U(  , ) 

2090  FOR  Jr2  TO  N-1 


1 
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SPLINFIT  (continued) 

2100  LET  N( J , 1 ) =L(  J) 

2110  LET  N( J,2)=2*(L( J)+L( J-1 ) ) 

2120  LET  N(J,3)=L(J-1) 

2130  FOR  K=1  TO  S 

2140  LET  B(K, J)=3*(L( J-1 )“2*(P(K, J+1 )-P(K, J) )+L( J)“2*(P(K, J)-P(K, J-1 ) ) ) 
2150  LET  B(K,J)=B(K,J)/(L(J)»L(J-1)) 

2160  NEXT  K 
2170  NEXT  J 

2180'  THE  FOLLOWING  IS  THE  GAUSSIAN  ELIMINATION  ON  THE  MATRIX, 

( 2190  FOR  1=2  TO  N 

1 2200  IF  N(I, 1 )=0  THEN  2310 

2210  LET  D=N(I-1 ,2)/N(I,1 ) 

2220  FOR  K=1  TO  3 

2230  LET  N(I,K)=N(I,K)*D-N(I-1 ,K+1 ) 

2240  LET  B(K,I)=B(K,I)*D-B(K,I-1 ) 

2250  NEXT  K 
2260  LET  Q=N(I,2) 

2270  FOR  K=1  TO  3 

2280  LET  N(I,K)=N(I,K)/Q 

2290  LET  B(K,I)=B(K,I)/Q 

2300  NEXT  K 

2310  NEXT  I 

2320  FOR  K=1  TO  S 

2330  FOR  J=0  TO  N-1 

2340  LET  U(K,N-J)=(B(K,N-J)-N(N-J,3)*U(K,N+1-J))/N(N-J,2) 

2350  NEXT  J 
2360  NEXT  K 
2370  SUBEND 

2380'  THE  BASIC  TANGENT  VECTORS  ARE  NOW  AVAILABLE  HERE. 

2390'  THE  FOLLOWING  MODIFIES  THE  VECTOR  MAGNITUDES  TO  EQUAL  THE 
2400'  SHORTER  OF  THE  CHORD  LENGTHS  TO  ADJACENT  POINTS  FOR  EACH  POINT. 
2410  SUB  "TANMOD":S,N,U(  , ),D(  ) 

2420  FOR  J=2  TO  N-1 

2430  IF  S=3  THEN  2460 

2440  LET  U=SQR(U( 1 , J)“2+U(2, J)“2) 

2450  GO  TO  2470 

2460  LET  U=SQR(U(1,J)'‘2+U(2,J)“2+U(3,J)"2) 

2470  IF  D( J-1 )<D(J)  THEN  2520 
- 2480  FOR  K=1  TO  S 

2490  LET  U(K,J)=(U(K,J)*D(J))/U 

2500  NEXT  K 

2510  GO  TO  2550 

2520  FOR  K=1  TO  S 

2530  LET  U(K,J)=(U(K, J)*D( J-1 ) )/U 

2540  NEXT  K 

2550  NEXT  J 

2560  LET  J=1 

2570  IF  S=3  THEN  2600 

2580  LET  U=SQR(U(  1 , J)”“2+U(2,  J)~2) 

2590  GO  TO  2610 
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SPLIMFIT  (continued) 

2600  LET  UsSQR(U(1lJ)"2+U(2lJ)'*2+U(3,J)'*2) 

2610  IF  J=N  THEN  2640 
2620  LET  FrD(J) 

2630  00  TO  2660 
2640  LET  F=D( J-1 ) 

2660  FOR  K= 1 TO  S 

2660  LET  d(K,J)=(J(K,J)*F)/U 

2670  ME XT  K 

2680  IF  J=N  THEM  2710 

2600  LET  J=N 

2700  GO  TO  2570 

2710  SUBEND 

2720'  THE  FOLLOWING  SOLVES  FOR  THE  PARAMETRIC  COEFFICIENTS  AND 
2730'  FILLS  IN  INTERMEDIATE  POINTS  BETWEEN  THE  SPECIFIED  POINTS, 

2740*  STORING  THE  DATA  IN  AN  OUTPUT  FILE  IN  STANDARD  GRAPHICAL 
2750'  FORMAT 

2760  SUB  "CURGEN":S,N,P(  , ),L(  ),Z,U(  , ),#1 
2770  PRINT  "DO  YOU  WANT  A LIST  OF  POINT  RESULTS"; 

2730  INPUT  F$ 

2790  IF  F$=" NO"  THEN  2830 
2600  PRINT 

2610  PRINT  "POINT" ;rA3(7);"F'(X)";TAB(2C);"F'(Y)";TAB(27);"F'(Z)"; 

2820  PRINT  TAB(37);"F"(X)";TAB(50);"F"(Y)";TAB(62);"F'  ’ (Z)" 

2830  FOR  J=1  TO  N-1 
2840  FOR  K=1  TO  S 
2850  LET  F( 1 , K)=P(K, J) 

2360  LET  F(2,fO=J(X,J) 

2870  LET  F(3,K)s(3/L(J)"2)*(P(K,J+1)-P('<,J))-(  1/L( J) )*(U(.K,  J+1  )+2*!J(X,  J) ) 
2830  LET  F(4 ,K)=(-2/L( J)"3) *( P(K , J+1 )-P(K , J) )+( 1/L( J)"2)*( U(K, J+1 )+J(K, J) ) 

£§00  NEXT  K 

2000  IF  F$="UO"  THEN  2070 
2010  IF  3=3  THEN  2050 

2020  PRINT  J ; TAB (7 ) ; F ( 2 , 1 ) ; TA3( 20 ) ; F ( 2 , 2 ) ; TAB( 27 ) ; 0 ; T AB( 37 ) ; 

2030  PRINT  2*F( 3, 1 ) ;TA3(50) ;2*F( 3,2) ;TAB(62) ;0 
20 4 j 10  TO  2070 

2050  PRINT  J;TAS(7);F(2,1);TAU(20);F(2,2);TAB(27);F(2,3);TAB(37); 

2060  PRINT  2*F(3,1);TA3(50);2*F(3,2);TAB(62);2*F(3,3) 

2070  FOR  T=0  TO  L(J)  STEP  L(J)/Z 

2030  IF  J=1  THEN  3010 

2000  IF  TOO  THEM  3010 

3000  GO  TO  3030 

3010  FOR  K=1  TO  S 

3020  LET  C(K)=F(l,iO+F(2,K)*T+F(3,K)»(T"2)+F(4,K)*(r3) 

3030  NEXT  K 

3040  IF  3=3  THEN  3070 

305°  PRINT  #1 :C( 1 ) ,";C(2) 

3060  GO  TO  3080 

'0  PRINT  :C( 1) ,";C(2) ,";C(3) 

>0  to  NEXT  T 
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SPLINFIT  (continued) 


IF  F$="N0"  THEN  3220 

PRINT  J+1 ;TAB(7) ;F(2, 1 )+2*F(3, 1 )*L(J)+3*F(4, 1 )*L( J)“2;TAB(20) ; 

PRINT  F(2,2)+2*F(3,2)*L(J)+3*F(4,2)*L(Jr2;TAB(27); 

IF  S=3  THEN  3160 
PRINT  0;TAB(37); 

GO  TO  3170 

PRINT  F(2,3)+2*F(3,3)*L(J)+3*F(4,3)*L(J)"2;TAB(37); 

PRINT  2«F(3,1)+6*F(4,1)*L(J);TAB(50);2»F(3,2)+6*F(4,1)»L(J);TAB(62); 
IF  S=3  THEN  3210 
PRINT  0 
GO  TO  3220 

PRINT  2*F(3,3)+6*F(4,3)*L(J) 

IF  S=3  THEN  3250 
PRINT  #1 :"1E37 , 1E37" 

GO  TO  3260 

PRINT  #1 :"1E37,1E37,1E37" 

SUBEND 

SUB  "GRAPH" :N,#1 ,S,D$ 

RESET  #1 

IF  S=2  THEN  3310 
CALL  "TRANS" :#1 
DIM  X(500),Y(500) 

LET  1=0 
LET  1=1+1 

IF  END#1  THEN  3370 
INPUT  #1 :X(I) ,Y(I) 

GO  TO  3330 
LET  N=I-1 

PRINT  "IS  THE  OUTPUT  TO  BE  DISPLAYED  ON  TELETYPE, CRT," 

PRINT  "TYPAGRAPH, COMPUTER, OR  LATER"; 

INPUT  Z$ 

IF  Z$="TELETYPE"  THEN  3480 

IF  Z$="CRT"  THEN  3500 

IF  Z$="TYPAGRAPH"  THEN  3520 

IF  Z$="COMPUTEK"  THEN  3540 

IF  Z$="LATER"  THEN  3560 

PRINT  "INCORRECT  FORMAT,  TRY  AGAIN" 

GO  TO  3380 

CALL  "TELSUB" : # 1 , 0 , A$ , 3 
GO  TO  3570 

CALL  "CRTSUB":#1,0,A$,3 
GO  TO  3570 

CALL  "TYPSUB" : # 1 , 0 , A$ , 3 
GO  TO  3570 

CALL  "COMSUB" : # 1 , 0 , A$ , 3 
GO  TO  3570 

PRINT  "POINTS  ARE  LOCATED  AS  SPECIFIED  IN  FILE  ";D$ 

SUBEND 

SUB  "TRANS" : 0 1 

PRINT  "INPUT  NAME  OF  SCRATCH  FILE"; 
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SPLINFIT  (continued) 

3600  INPUT  B$ 

3610  FILE  #2:B$ 

3620  SCRATCH  #2 

3630  PRINT  "TYPE  OF  PROJECTION: ( 1=X-Y,2rY-Z,3=Z-X)" 

3640  INPUT  A 

3650  ON  A GO  TO  3660,3710,3760 

3660  IF  MORE  #1  THEN  3680 

3670  GO  TO  3810 

3680  INPUT  #1 :A,B,C 

3690  PRINT  #2:A;",";B 

3700  GO  TO  3660 

3710  IF  MORE  #1  THEN  373 0 

3720  GO  TO  3810 

3730  INPUT  #1 :A,B,C 

3740  PRINT  #2:B;",";C 

3750  GO  TO  3710 

3760  IF  MORE  #1  THEN  3780 

3770  GO  TO  3810 

3780  INPUT  #1 :A,B,C 

3790  PRINT  #2:C;",";A 

3800  GO  TO  3760 

3810  SCRATCH  #1 

3820  RESET  #2 

3830  INPUT  #2:A,B 

3840  PRINT  #1:A;",";B 

3850  IF  END  #2  THEN  3870 

3860  GO  TO  3830 

3870  RESET  #1 

3880  SUBEND 

3890  SUB  "MATINV":S,N,P(  , ),N(  , ),B(  , ),D(  ),L(  ),U(  , ),A$,B$,Z 
3900  DIM  M(70,70),C(70),V(70,70),W(70) 

3910  IF  SEG$(A$,1,3)="CYC"  THEN  3980 
3920  IF  SEG$(A$,1,3)="ANT"  THEN  3980 
3930  MAT  M=ZER(N,N) 

3940  MAT  C=ZER(N) 

3950  MAT  W=ZER(N) 

3960  MAT  V=ZER(N,N) 

3970  GO  TO  4020 

3980  MAT  M=ZER(N-1 ,N-1 ) 

3990  MAT  C=ZER(N-1) 

4000  MAT  W=ZER(N-1 ) 

4010  MAT  V=ZER(N-1,N-1) 

4020  IF  SEG$(A$,1,3)="CYC"  THEN  4190 
4030  IF  SEG$(A$,1,3)="ANT"  THEN  4190 
4040  IF  SEG$(A$, 1 , 1 )="Q"  THEN  4080 
4050  LET  M( 1 , 1 )=N( 1 ,2) 

4060  LET  M(1,2)=N(1,3) 

4070  GO  TO  4110 
4080  LET  M( 1 , 1 )=N(1 , 1 ) 

4090  LET  M(1,2)aN(1,2) 
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SPLINFIT  (continued) 


4600  1AT  V=IiJV(  5) 

4510  If  3EG$(  A1>,  1 ,3)="CYC"  THEN  4730 
4620  IF  SEG$ ( A$ » 1 , 3) =" ANT"  THEtsJ  4730 
4530  FOR  K=1  TO  S 
4640  FOR  J=1  TO  .4 
4550  LET  C(J)=B(K,J) 

4660  NEXT  J 
4670  MAT  ,'/=V*C 
4680  FOR  J=1  TO  H 
4600  LET  U(K, J)  = ./( J) 

4700  MSXT  J 
4710  NEXT  K 
4720  GO  TO  4000 
4730  FOR  K=1  TO  3 
4740  FOR  J=1  TO  ii-1 
4750  LET  C( J)=8(K, J) 

4760  NEXT  J 
4770  MAT  W=V*C 
4730  FOR  J=1  TO  0-1 
4700  LET  LJ(K, J)=V( J) 

4300  NEXT  J 
4310  NEXT  K 

4320  IF  3EG$(A7>,  1 ,3)="ANT"  THEN  4370 
4330  FOR  K=1  TO  S 
4340  LET  U(K,N)=U(K, 1 ) 

4350  NFXT  !< 

4360  GO  TO  4000 
4370  FOR  K= 1 TO  S 
4830  LET  U(K,N)=-U(fC,1) 

4300  JEXT  X 
4000  SUBEND 
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COEFFICIENTS  FOR  EXPANSION  OF  BUBBLE-CASE  1 


AO 

A2 

A4 

A6 

1 « 

0.695397 

0.412924 

-0.140418 

5.00495  E-2 

l t 

0.632829 

0.355654 

-0.121983 

5.96184  E-2 

0.567656 

0.400168 

-0.123263 

3.51356  E-2 

1 - 

0.563461 

0.321668 

-9.88701  E-2 

4.54331  E-2 

s*  f 

f 1 

0.5126 

0.242358 

-5.75205  E-2 

3.68124  E-2 

0.487923 

0.225188 

-9.37463  E-2 

5.5635  E-3 

j 

0.509801 

0.230113 

-0.15733 

3.15432  E-2 

i 1 i 

0.536802 

9.79686  E-2 

-0.179822 

-1.19548  E-2 

r 

0.507123 

-5.00193  E-2 

-0.221311 

1.97178  E-2 

0.500335 

-8.56244  E-2 

-0.258396 

1.43054  E-3 

1 


M 


J 


COEFFICIENTS  FOR  EXPANSION  OF  BUBBLE-CASE  2 


AO 

0.681519 

A2 

0.444012 

A4 

-0.210683 

A6 

9.44173  E-2 

0.685496 

0.412312 

-0.210071 

0.096682 

0.68136 

0.354767 

-0.143659 

7.10024  E-2 

0.66308 

0.359718 

-0..  139868 

6.56617  E-2 

0.638041 

0.364293 

-0.111248 

5.17293  E-2 

0.590996 

0.426635 

-0.108948 

3.94992  E-2 

0.545938 

0.38489 

-0.091815 

5.19594  E-2 

0.520715 

0.313837 

-1.99842  E-2 

2.28483  E-2 

0.464311 

0.27543 

1.74285  E-2 

2.29671  E-2 

0.427852 

0.219103 

3.42362  E-2 

5.23286  E-2 

0.370364 

0.186245 

1.53092  E-3 

0.092001 

0.342797 

0.116947 

-3.43396  E-2 

6.32129  E-2 

0.33547 

0.040418 

-0.113346 

-3.14466  E-2 

0.359263 

-5.14157  E-3 

-0.139174 

-8.07007  E-2 

! 


CASF 

1 

CASF 

2 

FRAME 

x~? 

FRAME 

X~2 

1 

1*13169  F-3 

1 

1.51453  F-2 

o 

1.66658  F-3 

2 

0,02603 

3 

1.14931  E-4 

3 

2.23721  F-2 

4 

2.34727  F-4 

4 

1.30353  F-2 

r. 

4.78711  F-4 

5 

8.17185  E-3 

0 

9.4867  E-5 

6 

6.82436  E-7 

7 

5.74483  E-5 

7 

3.09257  E-2 

R 

5.69113  E-4 

8 

1.53104  E-2 

9 

1.96173  E-3 

9 

6.51845  E-3 

10 

1.45965  E-3 

10 

4.46423  E-3 

11 

7.06414  E-3 

1? 

6.78229  E-3 

13 

0.100918 

14 

0.16202 
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APPENDIX  A7 


Program  written  for  the  USNA  PDP-11,  which 
generates  a three  dimensional  display  for 
the  E+S  Picture  System. 
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100  OPEN  "AKCOEF.DAT"  FOR  INPUT  AS  FILE  #3 

110  INPUT  #3:M 

120  DIM  C(72,4) ,D(599) ,A(4) 

130  FOR  1=1  TO  M 

140  INPUT  #3:C(I,1),C(II2),C(I,3),C(I,4) 

150  NEXT  I 

160  CALL  "PSNT"(2,0) 

170  H= 16300 
180  LET  V=30 
190  V= 15 

200  V9= 1 .00000E-03 
210  LET  S9=. 125663 
220  S9=S9/2 
230  P=2#H 
240  LET  S=10 

250  CALL  "WIND" ( -H , H , -H , H , H , -H , P ) 

260  T=-1000 

270  CALL  "TRAN"(0,0,T) 

280  B=0 

290  PRINT  B 

300  CALL  "ROT"(B, 1 ) 

310  FOR  C=1  TO  M 
320  FOR  1=1  TO  4 
330  LET  A(I)=C(C,I) 

340  NEXT  I 

350  GOSUB  490 

360  FOR  A=0  TO  360  STEP  S 

370  CALL  "PUSH" 

380  CALL  "ROT" (A, 2) 

390  CALL  "DR3D"(D,INT(I/3),2,2) 

400  CALL  "POP" 

410  NEXT  A 
420  CALL  "NFRM" 

430  CALL  "SWCH"(0,Q9) 

440  IF  Q9<>1  THEN  470 
450  CALL  "RSTP" 

460  STOP 
470  NEXT  C 
480  STOP 
490  REM 
500  LET  1=0 

510  FOR  Q=0  TO  1.5708  STEP  S9 
520  LET  P0=1 

530  LET  P2=( 3*COS(2*Q)+1 )/4 

540  LET  P4=(35*C0S(4*Q)+20*C0S(2*Q)+9)/64 

550  LET  P6=(231#COS(6»0)+126»COS(4*Q)+105*COS(2«Q)+50)/512 

560  LET  R=A( 1 )*P0+A(2)*P2+A(3)*P4+A(4)*P6 

570  LET  X=R*SIN(Q) 

580  LET  Y=R*COS(Q) 

590  IF  Y<0  THEN  670 


APPENDIX  A8 

J 

Program  for  hard  copy  generation  of  three 
dimensional  display  from  PDP-ll/E+S  Picture 
System  on  Xynetics  plotter. 
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AKFTBD 

100  OPEN  "AKCOEF.DAT"  FOR  INPUT  AS  FILE  #3 

110  INPUT  #3:M 

120  DIM  C(72,4),D(599),A(4) 

130  REM 

140  FOR  1=1  TO  M 

150  INPUT  #3:C(I,1),C(I,2),C(I,3),C(I,4) 

160  NEXT  I 

170  CALL  "PSNT"(2,0) 

180  LET  H= 16000 
190  LET  V=1 
200  LET  V9=.2 
210  LET  S9=. 125663 
220  S9=S9/2 
230  LET  V8=0 
240  P=2*H 
250  LET  S=10 
260  LET  B=50 
270  FOR  C=1  TO  M 
280  FOR  1=1  TO  4 
290  LET  A(I)=C(C,I) 

300  NEXT  I 

310  CALL  " INST" ( -2000 , 2000 , 22000- V8 , 26000- V8 , 2000 , -2000 ) 
320  LET  V8=V8+5000 
330  GOSUB  520 

340  CALL  "WIND"(-H,H,-H,H,H,-H) 

350  CALL  "MAST"(-H,H,-H,H,H,-H) 

360  CALL  "ROT"(B, 1 ) 

370  FOR  A=0  TO  360  STEP  S 
380  CALL  "PUSH" 

390  CALL  "R0T"(A,2) 

400  CALL  "DR3D"(D,INT(I/3),2,2) 

410  CALL  "POP" 

420  NEXT  A 
430  CALL  "POP" 

440  NEXT  C 
450  CALL  "NFRM" 

460  CALL  "LITE"(0,1) 

470  CALL  "SWCH"(0,Q9) 

480  IF  Q9<>1  THEN  470 
490  CALL  "LITE" (0,0) 

500  CALL  "RSTP" 

510  STOP 
520  REM 
530  LET  1=0 

540  FOR  Q=0  TO  1.5708  STEP  S9 
550  LET  C2=C0S(2*Q) 

560  LET  C4=COS(4*Q) 

570  LET  C6=COS(6*Q) 

580  LET  P0=1 

590  LET  P2=(3#C2+1)/4 


-2- 

AKFTBD  (continued) 

600  LET  P4=(35*C4+20»C2+9)/64 
610  LET  P6 = ( 23 1 *C6+ 1 26*C4+ 1 05*C2+50 ) /5 1 2 
620  LET  R=A( 1 )#P0+A(2)*P2+A(3)*P4+A(4)*P6 
630  LET  X=R»COS(Q) 

640  LET  Y=R*SIN(Q) 

650  LET  Y=Y-.5 
660  IF  Y<0  THEN  740 
670  IF  ABS(X)<V9  THEN  740 
680  LET  D(I)=INT(X*1000*V) 

690  LET  1=1+1 

700  LET  D( I ) =INT ( Y* 1 000*V ) 

710  LET  1=1+1 
720  LET  D(I)=INT(0) 

730  LET  1=1+1 
740  NEXT  Q 
750  RETURN 
760  END 
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The  report  concerns  the  behavior  of  cyclic  liquid  jets  in  pulsating 


The  jets  were  developed  in  bubbles  which  were  situated  on  a vibrating 
platfom  that  was  contained  within  a vessel  capable  of  sustaining  reduced  pres- 
sures. The  entire  assembly  was  mounted  on  a variable  frequency  vibration  table 
and  driven  at  displacement  amplitudes  on  the  order  of  a millimeter.  The  reduced 
pressure  enabled  the  jets  to  be  quite  easily  formed  within  the  vessel  ; the  jets 
were  formed  in  bubbles  containing  both  vapor  and  gas,  thus  resulting  in  bubbles 
exhibiting  Jet  behavior  over  several  cycles.  The  bubbles  were  photographed  with 
a high  speed  camera  and  a detailed/analysis  made  of  the  behavior  of  the  Jet  as  a 
function  of  time. 

It  was  found  that  the  bubble  surface  could  be  well  represented  by  the 
first  four  even  terms  of  a Legendre  polyasmial  expansion.  Further,  a graphic 
three-dimensional  display  of  the  bubble  surface ^ generated  by  the  expans ion, gives 
excellent  agreement  with  the  observed  behavior.  3 

Data  of  calculations  is  given  in  the  appendices. 


